Four common vitamin D receptor polymorphisms and coronary artery disease susceptibility: A trial sequential analysis

Background Studies on the susceptibility of vitamin D receptor (VDR) polymorphisms to coronary artery disease (CAD) reached controversial results. We performed this study for a more accurate evaluation between the VDR polymorphisms and CAD susceptibility. Methods PubMed, Embase, CNKI, Wan Fang, and VIP databases were searched. The odds ratios (ORs) and 95% confidence intervals (95% CIs) were used to evaluate the associations. Trial sequential analysis (TSA) was introduced to estimate the positive associations. The potential functions of the VDR polymorphisms were analyzed based on the SNPinfo and ENSEMBL databases. Results Thirteen studies were finally included. In the overall analysis, increased CAD risks were observed in the VDR rs1544410 polymorphism and verified by the TSA; for the rs2228570 and rs731236 polymorphisms, significant associations with high heterogeneity were detected; decreased risk was remarkably observed for the rs7975232 polymorphism. In the subgroup analysis, wide associations with reduced heterogeneity were observed in the rs2228570, rs1544410, and rs731236 polymorphisms. The RNAfold analysis indicated the mutant G allele of the rs1544410 polymorphism was easier to disperse from the DNA double helix structure and may have a potential crucial role in the VDR transcription process. Conclusions Our analysis supports the role of the rs1544410 polymorphism in the VDR gene as a risk factor for CAD. The VDR rs2228570 and rs731236 polymorphisms were associated with increased CAD risks in the White population. Restrict decreased CAD risk was firstly discovered in the rs7975232 polymorphism. Limitations Firstly, the language was restricted to English and Chinese, which will cause the limited number of studies included; secondly, other unknown polymorphisms in VDR polymorphisms could also be associated the CAD susceptibility, and more case-control studies with comprehensive clinical outcomes and GWAS studies were required; thirdly, the rs1544410, rs7975232 and rs731236 polymorphism are in strong LD, haploid factors with CAD risk need to be considered; fourthly, the mechanisms of the VDR polymorphism on the VDR gene or RNA or protein were not discussed enough, further mechanistic studies are required; at last, genetic factor was the one side for CAD susceptibility, the interaction between environmental risk factors should be considered.


Introduction
Risk factors have always been a hot topic in the study of coronary artery disease (CAD), which is the main leading cause of death in the world [1][2][3]. Traditional risk factors, such as smoking, obesity, high blood lipids, etc. [4], help physicians guide the population to prevent CAD. Whether genetic factors have an influence on the CAD risk remains unclear. A recent study reported the heritability of CAD has been estimated between 40% and 60% [5], which implied that genetic factors would play a distinctive role in the susceptibility of CAD.
Reduced serum vitamin D concentration was reported to be an increased risk marker for CAD [6], Vitamin D receptor (VDR) is a vital signal transduction molecule for vitamin D [7]. The VDR gene is located on chromosome 12q13.1, and has four common single nucleotide polymorphisms (SNPs) which are rs2228570 (FokI F/f in exon 2), rs1544410 (BsmI B/b in intron 8), rs7975232 (ApaI A/a in intron 8) and rs731236 (TaqI T/t in exon 9) [8]. Van Schooten et al. firstly reported the rs1544410 polymorphism was associated with the severity of CAD [9]. A small group study then conducted by Ortlepp et al. also confirmed the former results in 2001 [10], but in a larger population study reported by him in 2003, no association was detected [11]. In the next decades, many studies were designed and conducted not only in the rs1544410 polymorphism but also in other VDR polymorphisms, however the conclusions were inconsistent.
We considered the inconsistence may owe to the bias in sample size, different characteristics of research population or the unavoidable system errors, therefore, a comprehensive study based on rigorous inclusion and exclusion criteria was performed, and trial sequential analysis was introduced to reduce the system errors and confirm our positive results, moreover, the function of the VDR polymorphisms were analyzed with bioinformatic tools. searching strategy. Each author independently reviewed the potential studies and the divergence were discussed by group-meeting held by Hailong Tao (The corresponding author).

Inclusion and exclusion criteria
Studies included in our study must meet the following inclusion criteria: (1) evaluation of the associations between the VDR polymorphisms and coronary artery disease susceptibility; (2) case-control study or cohort design; (3) detailed genotype frequency data could be acquired to calculate odds ratios (ORs), 95% confidence intervals (CIs) and Hardy Weinberg Equilibrium test; Exclusion criteria: (1) duplication of previous publications; (2) comment, review, case reports, animal studies and editorials; (3) study with no detailed genotype frequency data. The first two authors conducted the selection of potential included studies independently based on the inclusion and exclusion criteria. Any disagreement was solved by a discussion with the corresponding author.

Data extraction
For each study, the following data were independently extracted by the first two authors and the corresponding author used a standardized form: first author's last name, year of publication, study country, region, age, BMI, 25(OH) vitamin D, genotyping methods, detail genotype frequency data of cases and controls, genotype distribution in CAD (coronary artery disease) populations and controls, quality score and the result of Hardy-Weinberg Equilibrium test.

Quality score assessment
The modified Newcastle-Ottawa scale (NOS) was used to evaluate the quality of included studies in our study (S2 Table) [13][14][15][16]. Each included study was scored and regarded as either low quality (score � 6) or high quality (score > 6) based on items such as the definition of representativeness of cases, source of controls, sample size, quality control of genotyping method, and Hardy-Weinberg equilibrium.

Statistics analysis
Review Manager, Version 5.3 (The Nordic Cochrane Centre, The Cochrane Collaboration; Copenhagen, Denmark) and STATA 12.0 (STATA Corp, LP) were used for all analysis. P < 0.05 was considered to be significant. Hardy-Weinberg equilibrium (HWE) was evaluated for each study by Chi-square test in control groups, and P < 0.05 was considered as a significant departure from HWE. Odds ratio (OR) and 95% confidence intervals (CIs) were calculated. The odds ratio (OR) is the ratio of odds of an event in one group versus the odds of the event in the other group. An odds ratio (OR) of 1.0 indicates that there is no difference in risk (or odds) between the groups being compared. An OR of more than 1.0 indicates an increase in risk (or odds) among the exposed compared to the unexposed, whereas an OR <1.0 indicates a decrease in risk (or odds) in the exposed group [17] Heterogeneity was evaluated by Q statistic (significance level of P < 0.1) and I 2 statistic (greater than 50% as evidence of significant inconsistency). If the P value for heterogeneity was >0.10 and I2 <50%, indicating an absence of heterogeneity between studies, the fixed-effects model (the Mantel-Hasenszel method) would be used; If the P value for heterogeneity was �0.10 or I2 �50%, indicating a high exist of heterogeneity between studies, and the randomeffects model (the DerSimonian and Laird method) would be used. Besides, subgroup analyses were stratified by Race (White, Asian and African), Hardy-Weinberg equilibrium (in accordance with HWE, departure from HWE), sample size (�500, <500), genotyping method (PCR-RFLP, PCR-Taqman, PCR-ABD), and random-effects model were applied in subgroup analysis for more conservative results. We applied the Bonferroni method, which controls for the false discovery rate (FDR), to adjust for multiple comparisons.

Trial sequential analysis (TSA)
TSA (The Copenhagen Trial Unit, Center for Clinical Intervention Research, Denmark) is a methodology that combines an information size calculation (cumulated sample sizes of all included trials) for a meta-analysis with the threshold of statistical significance (http://www.ctu. dk/tsa). If the data in included studies are sparse or if there is repeated testing for significance in conducting an updated meta-analysis, the type I errors and type II errors are unavoidable [18,19].
To reduce the risk of type I errors, TSA was introduced in our analysis. The required information size was calculated according to an overall type-I error of 5%, a power of 80% and a relative risk reduction (RRR) assumption of 20% [20,21]. A continuity correction of 0.5 was also applied in zero-event trials [22].

Bioinformatics analysis
Ensembl is a genome browser for vertebrate genomes that supports research in comparative genomics, evolution, sequence variation and transcriptional regulation, and this database provides the genomic context, genes and regulatory elements, flanking sequence, population genetics, phenotype data, sample genotypes, linkage disequilibrium and phylogenetic context of a single nucleotide polymorphism (http://asia.ensembl.org/index.html). SNPinfo is an important bioinformatics analysis tool that predicts SNP function. The SNPinfo database can help researchers specify genes or linkage regions and select SNPs based on GWAS results, calculate linkage disequilibrium (LD), and predict functional characteristics of both coding and non-coding SNPs (https://snpinfo.niehs.nih.gov/) [23]. In addition, the RNAfold web server is one of the core programmes of the Vienna RNA package that has been used to predict the minimum free energy of single sequences that influence the stability of the structure [24]. Therefore, we conducted bioinformatics analyses using the above databases and methods to identify the potential molecular mechanisms for further research.

Results
The PRISMA flow diagram of the literature selection process was showed in S1 Table.

The analysis of VDR polymorphisms and CAD susceptibility
The results of overall and subgroup populations were showed in Table 2.      and African. As for the subgroup analyses based on HWE, Sample size and Genotyping method, wide significant associations with unsolved heterogeneity were observed.

rs7975232 polymorphism and CAD susceptibility
Interestingly, decreased risks of CAD were firstly discovered in overall analysis and subgroup analysis based on Race and Sample size. In overall analysis, decreased CAD risks were detected in allelic (OR [95%CI] = 0.93[0.88, 1.00], I2 = 4%) (Fig 1D)

Sensitivity analysis of associations between VDR polymorphisms and CAD susceptibility
We conducted the sensitive analyses on VDR polymorphism and CAD risk by omitting one study at a time in the calculation of the summary outcome (Fig 3). The results showed that no single study fundamentally changed the associations between these four VDR polymorphisms and CAD risk, which indicated that our meta-analysis results were relatively stable.

Publication bias
The Egger's test was introduced to analyze the publication bias, the P value for the test of these four VDR polymorphisms were 0.423 (rs2228570), 0.218 (rs1544410), 0.396 (rs731236) and 0.980 (rs7975232), respectively. Moreover, the Begg's funnel plots of these four polymorphisms were symmetrical (Fig 4). The results based on the Egger's test and the Begg's funnel plots indicated no publication bias for these four VDR polymorphisms with CAD risk.

Trial sequential analysis of associations between VDR polymorphisms and CAD susceptibility
Based on our analysis, increased CAD risks with low heterogeneity in the overall analysis of the rs1544410 polymorphism and the White population of the rs2228570 and rs731236 polymorphisms were discovered. Therefore, a trial sequential analysis was introduced to validate that our discoveries above were not false positive results. The allelic genetic model is a natural model of inheritance with a stronger genotype-phenotype association, which also does not pre-assume any interactions between the numbers of variant alleles. Therefore, we chose the allelic genetic model of the rs1544410 polymorphism in overall population and the rs2228570 and rs731236 polymorphism in the White population to conduct the trial sequential analysis. The X and Y axes represent the number of patients and the cumulative Z score, respectively. Within the designed assumptions of confidence and effect size, the information size for the rs1544410 polymorphism are 152472, the Z curves not only cross the statistical significance line (Z = 1.96, P = 0.05), but also cross the O' Brien Fleming boundaries (Fig 5), indicating that the significance level of our study was a true positive result. However, for the rs2228570 and rs731236 polymorphisms, although the Z curves cross the statistical significance line (Z = 1.96, P = 0.05), but not cross the O' Brien Fleming boundaries, which indicated more studies were required, and the information size for the rs2228570 and rs731236 polymorphisms were 84534 and 42415 respectively (Fig 6).

Bioinformatics analysis
Based on the genomic context obtained from the Ensembl database, the VDR rs2228570 polymorphism caused a "start lost", the rs1544410 and rs7975232 polymorphisms were intron variants, the rs731236 polymorphism was the synonymous variant. Hence, we analyzed the sequences of the four polymorphisms and the results from the SNPinfo database showed the VDR rs2228570 and rs731236 polymorphism were predicted the function of Splicing (Table 3). In addition, the secondary structure of DNA at the VDR rs1544410 sequences was predicted using RNAfold. The minimum free energy (MFE) and the free energy of the thermodynamic ensemble (FETE) of the rs1544410 polymorphism were -264.30 kcal/mol and -276.99 kcal/mol for the wild A allele, -265.80 kcal/mol and -278.74 kcal/mol for the mutant G allele, respectively. Based on the predicted free energy of the rs1544410 polymorphisms, the secondary structure of the polymorphisms was determined. Compared to the wild allele, the mutant alleles of the rs1544410 polymorphism caused a structure change which was pointed with arrows in Fig 7.

Discussion
Coronary artery disease (CAD) is a disease with very high morbidity and mortality. Early prevention based on genetic polymorphism can reduce the incidence of CAD [37,38]. In our study, four common single nucleotide polymorphisms (SNPs) in vitamin D receptor (VDR) gene (rs2228570, rs1544410, rs731236 and rs7975232) were comprehensively analyzed and subgroups analysis based on race, samples size, genetic features were performed. There was no genome-wide association study regarding the associations between the VDR polymorphisms and CAD susceptibility, and the four common VDR polymorphisms were widely discussed   with inconsistent results, therefore, we chose these four common VDR polymorphisms to investigate the associations between the four VDR polymorphisms and coronary artery disease (CAD) susceptibility.
In the previous meta-analysis [39][40][41], increased risks were found both in rs1544410 and rs731236 polymorphism, which agreed with our results, but decreased risk of rs2228570 was observed in his study. After careful analysis, we supposed the small sample size and different data recruiting methods could contribute to the discrepancy. In the data retrieving of two included studies (Ferrarezil et al. [35] and Nezhad et al. [34]), the author pooled different small group based on severity of CAD into one group, which caused extremely high heterogeneity. In our analysis, we extracted each small group as one single study to reduce the heterogeneity and positive results with no or subtle heterogeneity were widely observed. Jiang L reported a dose-response meta-analysis based on full subgroups stratified by sex, age, race, et al. and found prospective evidence for further testing of the utility of ferritin levels in predicting T2D

PLOS ONE
The vitamin D receptor polymorphisms and coronary artery disease susceptibility risk in a sex-specific manner [42,43], therefore, we completed exhaustive subgroup analysis stratified by race, HWE, sample size and genotyping method to explore the source of heterogeneity and the potential associations in subgroup, and many interesting results were discovered.
In the overall analysis, the rs1544410 polymorphism was discovered to be associated with an increased risk of CAD in all five genetic models and the positive results were verified by TSA in the allelic genetic model, which indicate that the role of the rs1544410 polymorphism in the VDR gene as a risk factor for CAD. For the mutant b allele, it has a 15% increased CAD risk compared to the B allele. In terms of genotype, the Bb and bb genotype have 22% and 29% increased CAD risk compared to the BB genotype, respectively. As for the rs2228570 and rs731236 polymorphism, increased risk with high heterogeneity were widely observed, the mutant f allele has a 27% increased risk compared to the F allele for the rs2228570 polymorphism; for the rs731236 polymorphism, the mutant t allele has a 19% increased risk compared to the T allele. Interesting findings emerged on the rs7975232 polymorphism, decreased risks were firstly observed. The mutant a allele has a 7% decreased CAD risk compared to the A allele, and the aa genotype has a 13% decreased CAD risk compared to the AA genotype, however, the relative small sample size could have an influence on the evaluation of the rs7975232 polymorphism, and more well-designed studies were required to solidate the potential protective role of the rs7975232 polymorphism.
In the subgroup analysis, the high heterogeneity of the rs2228570 and rs731236 polymorphisms were significantly reduced in the White population. In the White population, increased CAD risks were extensively detected in the rs1544410, rs2228570 and rs731236 polymorphism. However, in Asian subgroup, an 60% increased CAD risk of the ff genotype is observed in the rs2228570 polymorphism compared to the Asian FF genotype, and no association is detected in the rs1544410 and rs731236 polymorphisms. The results in subgroup analysis stratified by race may indicate the White population with rs2228570, rs1544410 and rs731236 are more susceptible to CAD. As for the rs7975232 polymorphism, in the White population, the mutant a allele has a 7% decreased CAD risk compared to the A allele. Sample size is an important parameter in the case-control studies. In the subgroup analysis based on sample size, we detected that increased or decreased risks of the four VDR polymorphisms were widely observed in the large subgroup, which implied case-control studies with sufficient sample size could discover more meaningful data. Homogeneity is a crucial factor in the statistical Analysis, therefore HWE and genotyping typing method were analyzed in different subgroups, however, the results showed these two-subgroup analysis did not seem to affect the high heterogeneity. Compared to the traditional risk factors like smoking, being overweight, and lack of exercise et al., the VDR polymorphisms associated with CAD susceptibility we found could help the population identify CAD earlier and provide individualized treatment.
Lower plasma level of vitamin D was associated with increased risk of CAD [6,[44][45][46]. VDR is the crucial signal transduction molecule in the vitamin D pathway. From animal research reported by Xiang et al. [47], overexpressing the vitamin D receptor could inhibit the formation of atherosclerotic plague in APOE-deficient mice. In the CAD population, the TT genotype of rs2228570 polymorphism had a lower serum level of vitamin D compared to CC genotype [29]; for rs1544410 and rs731236 polymorphisms, the mutant genotype was associated with the lower plasma level of vitamin D [34,35]. The polymorphisms in VDR may have an influence in the interaction between VDR and Vitamin D and the serum level of Vitamin D. Causal inference analysis analyze the functional polymorphisms in a gene whether can causally trigger the development of a related disease through mediating the expression of this gene in specific tissues [48,49]. Zhang F et al. reported the genetically determined PTSD confers a causal effect on depression and depressed affect, but not major depressive disorder [50], moreover, deep learning or machine learning is a hot topic in classification and prediction of diseases based on biomarkers [51,52], which inspired us to conduct the causal inference analysis of the functional VDR polymorphisms in CAD and discuss the possibility to use the vitamin D receptor genetic variants related to CHD for the prediction or early diagnosis of CHD in our next mechanism study.
The VDR rs2228570 polymorphism caused a "start lost", the rs1544410 and rs7975232 polymorphisms were intron variant, the rs731236 polymorphism was the synonymous variant. The VDR rs2228570 and rs731236 polymorphism were predicted the function of Splicing. In addition, the secondary structure of DNA at the VDR rs1544410 sequences was predicted by using the RNAfold, which indicated that the mutant G allele could cause an easier dispensation from the DNA double helix structure. The SNP in rs2228570 polymorphism is located in the exon 2, which is near the translation start sequence, and the mutant T allele causes a structural modification of three amino acids longer protein leading to the change of potential protein function [53]. Unlike the rs2228570 polymorphism, the rs1544410, rs7975232 and rs731236 polymorphisms are located near the 3' end of the gene and cause no structural transformation [53], but they have a strong linkage disequilibrium (LD) [54]. The AAC haplotype composed by the A allele of rs1544410, A allele of rs7975232 and C allele of rs731236 was associated with an increased risk of CAD in type 2 diabetes subjects reported by Ferrarezi et al. [35], furthermore, a VDR GATG haplotype (G allele of rs731236, A allele of rs7975232, T allele of rs1544410 and G allele of rs2228570) was found to be associated with atherosclerotic disease in rheumatoid arthritis patients [55]. These studies suggest a joint role of the three polymorphisms in CAD susceptibility. Besides, acetyl-cytidine on RNA expression is also playing key role on the human diseases. Gehui Jin et al. reported the role and mechanism of ac4C in geneexpression regulation and demonstrated the relevance of ac4C to a variety of human diseases [56]. We found the changed RNA second structure in mutant allele of the VDR polymorphism, the changed structure may provide potential acetyl-cytidine loci and affect the RNA expression, which provide direction for our next mechanistic studies.
There were several limitations in our meta-analysis. Firstly, the language was restricted to English and Chinese, which will cause the limited number of studies included; secondly, other unknown polymorphisms in VDR polymorphisms could also be associated the CAD susceptibility, more case-control studies with comprehensive clinical outcomes and GWAS studies were required; thirdly, the rs1544410, rs7975232 and rs731236 polymorphism are in strong LD, haploid factors with CAD risk need to be considered; fourthly, the mechanisms of the VDR polymorphism on the VDR gene or RNA or protein were not discussed enough, further mechanistic studies are required; at last, genetic factor was the one side for CAD risk, the interaction between environmental risk factors should be considered.

Conclusion
Our analysis supports the role of the rs1544410 polymorphism in the VDR gene as a risk factor for CAD. The VDR rs2228570 and rs731236 polymorphisms were associated with increased CAD risks in the White population. Restrict decreased CAD risk was firstly discovered in the rs7975232 polymorphism.
Supporting information S1